###############################
# Overall ATT: Crime Concerns
###############################

rm(list=ls())

library(Hmisc)
library(readstata13)
library(foreign)
library(tidyverse)
library(did)
library(ggplot2)
set.seed(1111)

#sink("05_overall_att_concerns.txt")

############################
# Read and prepare data 
############################

# Load point estimates and standard errors
d1 <- read.table("overallATT_1.txt",comment.char="")
d2 <- read.table("overallATT_2.txt",comment.char="")
d3 <- read.table("overallATT_3.txt",comment.char="")
d4 <- read.table("overallATT_4.txt",comment.char="")

# Prepare multiple comparisons
con_all = d1$V1[1]/d1$V1[2]
con_all_cov = d2$V1[1]/d2$V1[2]
con_haiti = d3$V1[1]/d3$V1[2]
con_haiti_cov = d4$V1[1]/d4$V1[2]

pv_con_all = 2*pnorm(q=con_all, lower.tail=T)
pv_con_all_cov = 2*pnorm(q=con_all_cov, lower.tail=T)
pv_con_haiti = 2*pnorm(q=con_haiti, lower.tail=F)
pv_con_haiti_cov = 2*pnorm(q=con_haiti_cov, lower.tail=F)

p1 <- c(pv_con_all,pv_con_all_cov)
p2 <- c(pv_con_haiti,pv_con_haiti_cov)

# Corrected p-values visas from all countries on concerns (Table A17)
round(p.adjust(p1, "bonferroni"),3)

# Corrected p-values visas from Haiti on concerns (Table A17)
round(p.adjust(p2, "bonferroni"),3)

# Prepare Table 17_concerns
pe_haiti = c(d3[1,], d4[1,])
pe_all = c(d1[1,], d2[1,])
se_haiti = c(d3[2,], d4[2,])
se_all = c(d1[2,], d2[2,])

round(pe_haiti,3)
round(se_haiti,3)
round(pe_all,3)
round(se_all,3)

# Table A17 concerns
results_concerns <- rbind(round(pe_haiti,3),round(se_haiti,3),round(pe_all,3),round(se_all,3)) 
column <- rbind("visas from haiti","","visas from all countries","") 
results_concerns <- cbind(column,results_concerns)
write.table(results_concerns, "tableA17_concerns.txt", sep="\t")

sink()
